1 

Stochastic Surveillance Strategies 
for Spatial Quickest Detection 

Vaibhav Srivastava Fabio Pasqualetti Francesco Bullo 



Abstract — We design persistent surveillance strategies for the 
quickest detection of anomalies taking place in an environment 
of interest. From a set of predefined regions in the environment, 
a team of autonomous vehicles collects noisy observations, which 
a control center processes. The overall objective is to minimize 
detection delay while maintaining the false alarm rate below 
a desired threshold. We present joint (i) anomaly detection 
algorithms for the control center and (ii) vehicle routing policies. 
For the control center, we propose parallel cumulative sum 
(CUSUM) algorithms (one for each region) to detect anomalies 
from noisy observations. For the vehicles, we propose a stochastic 
routing policy, in which the regions to be visited are chosen 
according to a probability vector. We study stationary routing 
policy (the probability vector is constant) as well as adaptive 
routing policies (the probability vector varies in time as a 
function of the likelihood of regional anomalies). In the context 
of stationary policies, we design a performance metric and 
minimize it to design an efficient stationary routing policy. Our 
adaptive policy improves upon the stationary counterpart by 
adaptively increasing the selection probability of regions with 
high likelihood of anomaly. Finally, we show the effectiveness of 
the proposed algorithms through numerical simulations and a 
persistent surveillance experiment. 

Index Terms — vehicle routing, statistical decision making, 
quickest detection, persistent surveillance, patrolling, security, 
motion planning. 

I. Introduction 

Recent years have witnessed a surge in the application of 
autonomous agents in various activities such as surveillance 
and information collection. In view of the recent Icelandic ash 
problem, the oil spill in the gulf of Mexico, and recurring wild 
fires, surveillance strategies resulting in the quickest detection 
of anomalies are of considerable importance. Due to extreme 
sensor and modeling uncertainties in these situations, robust 
anomaly detection methods need to be employed. Generally, 
a limited number of vehicles are deployed to survey a large 
number of regions, and it is fundamental that the vehicles 
collect the information that is most effective to minimize the 
detection delay of anomalies. In this paper we design surveil- 
lance strategies that result in quick detection of anomalies. 

A preliminary version of this work ( Srivastava and Bullo 2011 1 was 
presented at IEEE Conference on Decision and Control and European Control 
Conference, 2011. In addition to the ideas in ( Srivastava and Bullo 2011), 
this paper contains a rigorous analysis of the single vehicle surveillance, 
the multiple vehicle surveillance, extensive numerical illustrations, and a 
persistent surveillance experiment. 
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07-1-0528, by NSF Award CPS-1035917 and by ARO Award W911NF-11- 
1-0092. 
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A reliable detection of anomalies can be achieved by col- 
lecting observations sequentially until the evidence suggesting 
an anomaly reaches a substantial level. Various sequential 
algorithms for the detection of anomalies have been presented 
in (Basseville and Nikiforov| [1993 ). Furthermore, it is known 
that a human being typically performs well in detecting and 
identifying anomalies from observations. Recent advances in 
cognitive psychology ( Bogacz et al. 2QQ6[ Ratnam et al.| 
|2QQ3| ), show that human performance in decision making tasks 
is well modeled by sequential statistical procedures such as 
the CUSUM algorithm. Inspired by the above human decision 
making models, in this work we adopt sequential statistical 
tests for anomaly detection. 

Our setup and approach. We consider an environment 
comprising of potentially disjoint regions of interest, and we 
employ a team of autonomous vehicles for the persistent 
surveillance of these regions. In particular, the vehicles visit 
the regions, collect information, and send it to a control center. 
We study a spatial quickest detection problem with multiple 
vehicles, that is, the simultaneous quickest detection of anoma- 
lies at spatially distributed regions when the observations for 
anomaly detection are collected by autonomous vehicles. For 
this problem, we let the control center run parallel CUSUM al- 
gorithms (one for each region) with the collected information. 
The control center then decides on the presence of anomalies 
in the regions. Finally, we design vehicle routing policies to 
collect observations at different regions. Our vehicle routing 
policies aim to minimize the anomaly detection time at the 
control center. 

Related work. Vehicle routing policies have witnessed a lot 
of attention in the robotics and controls literature. A survey 
on dynamic vehicle routing policies for servicing tasks is 
presented in ( [Bullo et al.[ |2Q11| ). Recently, the routing for 
information aggregation has been of particular interest. Klein| 
l et al.| (|2010') present a vehicle routing policy for optimal 
localization of an acoustic source. They consider a set of 
spatially distributed sensors and optimize the trade-off between 
the travel time required to collect a sensor observation and the 
information contained in the observation. They characterize 
the information in an observation by the volume of the Cramer- 
Rao ellipsoid associated with an optimal estimator. |Hollinger| 
|et al.| ( |2Qllb| ) study routing for an AUV to collect data from 
an underwater sensor network. They developed approximation 
algorithms for variants of the traveling salesperson problem 
to determine efficient policies that maximize the information 
collected while minimizing the travel time. Gupta et al.| ( |2006| ) 
study the estimation in a linear dynamical system with the 
observations collected by a set of mobile sensors. They deter- 
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mine stochastic trajectories for mobile sensors that minimize 
the error covariance of the Kalman filter estimate. Zhang et al.' 
( 2011| ) study the estimation of environmental plumes with 
mobile sensors. They minimize the uncertainty of the estimate 
of the ensemble Kalman filter to determine optimal trajectories 
for a swarm of mobile sensors. 

There has been some interest in decision theoretic infor- 
mation aggregation and vehicle routing as well. ICastanon' 
( fl995 ) poses the search problem as a dynamic hypothesis 
test, and determines the optimal routing policy that maximizes 
the probability of detection of a target. |Chung and Burdick 



( 2Q12| ) study the probabilistic search problem in a decision 
theoretic framework. They minimize the search decision time 
in a Bayesian setting. Certain optimal information aggrega- 
tion strategies for sequential hypothesis testing have been 
developed in ( [Srivastava et al.[ |2Qlla|b| ). [HoUinger et al. 



(|2011a) study an active classification problem in which an 
autonomous vehicle classifies an object based on multiple 
views. They formulate the problem in an active Bayesian 
learning framework and apply it to underwater detection. The 
persistent surveillance problem in this paper also concerns 
with decision-theoretic information aggregation and vehicle 
routing. In contrast to the aforementioned works that focus 
on classification or search problems, our focus is on quickest 
detection of anomalies. 

The problem of surveillance has received considerable atten- 
tion recently. Preliminary results on this topic have been pre- 



sented infChevaleyre , 2004 i Elmalia ch et al.[ |2008[ [Kingston 
|et al.| |2008| ). jPasqualetti et al.| ( |2012b| ) study the problem of 
optimal cooperative surveillance with multiple agents. They 
optimize the time gap between any two visits to the same 
region, and the time necessary to inform every agent about 
an event occurred in the environment. [Smith and Rus| ([2010|) 
consider the surveillance of multiple regions with changing 
features and determine policies that minimize the maximum 
change in features between the observations. A persistent 
monitoring task where vehicles move on a given closed path 
has been considered in ( jPasqualetti et al.j |2012a[ [Smith et aL] 
2Q12[ ), and a speed controller has been designed to minimize 
the time lag between visits of regions. 

Stochastic surveillance and pursuit-evasion problems have 
also fetched significant attention. In an earlier work, |Hespanha[ 



et al.j ( [1999| ) studied multi-agent probabilistic pursuit evasion 
game with the policy that, at each instant, directs pursuers 
to a location that maximizes the probability of finding an 
evader at that instant. [Grace and Baillieul| ( |2005| ) formulate 
the surveillance problem as a random walk on a hypergraph 
and parametrically vary the local transition probabilities over 
time in order to achieve an accelerated convergence to a 
desired steady state distribution. [Sak et aL] 02008) present 



partitioning and routing strategies for surveillance of regions 



for different intruder models. jSrivastava et al.j ( [2009] ) present 
a stochastic surveillance problem in centralized and decen- 
tralized frameworks. They use a Markov chain Monte Carlo 
method and a message passing based auction algorithm to 
achieve the desired surveillance criterion. They also show 
that the deterministic strategies fail to satisfy the surveillance 
criterion under general conditions. In this paper, we focus on 



stochastic surveillance policies. In contrast to aforementioned 
works on stochastic surveillance that assume a surveillance 
criterion is known, this work concerns the design of the 
surveillance criterion. The policies designed in this paper 
direct a vehicle with high probability to a region with high 
probability of being anomalous, a feature akin to the heuristic 
policy for the pursuer in ([Hespanha et al.j 1999] ). On the 
other hand, with respect to (Hespanha et al."'1999), our policy 
takes into account environmental factors, e.g., travel times and 
detection difficulty, and it satisfies an optimality criterion. 

Paper contributions. The main contributions of this work 
are fivefold. First, we formulate the stochastic surveillance 
problem for spatial quickest detection of anomalies (Section 
[n]). We propose the ensemble CUSUM algorithm for a control 
center to detect concurrent anomalies at different regions from 
collected observations (Section |III|. For the ensemble CUSUM 
algorithm we characterize lower bounds for the expected 
detection delay and for the average (expected) detection delay 
at each region. Our bounds take into account the processing 
times for collecting observations, the prior probability of 
anomalies at each region, and the anomaly detection difficulty 
at each region. 

Second, for the case of stationary routing policies, we pro- 
vide bounds on the expected delay in detection of anomalies 
at each region (Section [IV]). In particular, we take into account 
both the processing times for collecting observations and the 
travel times between regions. For the single vehicle case, we 
explicitly characterize the expected number of observations 
necessary to detect an anomaly at a region, and the corre- 
sponding expected detection delay. For the multiple vehicles 
case, we characterize lower bounds for the expected detection 
delay and the average detection delay at the regions. As a 
complementary result, we show that the expected detection 
delay for a single vehicle is, in general, a non-convex function. 
However, we provide probabilistic guarantees that it admits a 
unique global minimum. 

Third, we design stationary vehicle routing policies to 



collect observations from different regions (Section IV). For 
the single vehicle case, we design an efficient stationary policy 
by minimizing an upper bound for the average detection delay 
at the regions. For the multiple vehicles case, we first partition 
the regions among the vehicles, and then we let each vehicle 
survey the assigned regions by using the routing policy as 
in the single vehicle case. In both cases we characterize the 
performance of our policies in terms of expected detection 
delay and average (expected) detection delay. 

Fourth, we describe our adaptive ensemble CUSUM algo- 
rithm, in which the routing policy is adapted according to the 
learned likelihood of anomalies in the regions (Section [V]). We 
derive an analytic bound for the performance of our adaptive 
policy. Finally, our numerical results show that our adaptive 
policy outperforms the stationary counterpart. 

Fifth and finally, we report the results of extensive numerical 
simulations and a persistent surveillance experiment (Sections 
[Vl|and [VlI| ). Besides confirming our theoretical findings, these 
practical results show that our algorithm are robust against 
realistic noise models, and sensors and motion uncertainties. 
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Fig. 1. Persistent Surveillance Setup. A set of n regions is surveyed by 
m < n vehicles. Each vehicle visits the regions according to some policy 
and collects evidence from the visited region. The collected evidence is sent to 
an anomaly detection algorithm. The anomaly detection algorithm processes 
the collected evidence and decides on the presence of an anomaly. It also 
provides the likelihood of an anomaly being present, which in turn is used by 
the vehicle routing algorithm. The anomaly detection algorithm and vehicle 
routing algorithm constitute the control center, which can be implemented 
on-board of a vehicle. 



II. Problem Setup 

We consider the persistent surveillance of a set of n disjoint 
regions with a team of m < n identicaQ autonomous vehicles 
capable of sensing, communicating, and moving from one 
region to another. In persistent surveillance, the vehicles visit 
the regions according to some routing policy, collect evidence 
(sensor observation), and send it to a control center. The 
control center runs an anomaly detection algorithm with the 
evidence collected by the vehicles to determine the likelihood 
of an anomaly being present at some region (the control 
center declares an anomaly if substantial evidence is present). 
Finally, the control center utilizes the likelihood of an anomaly 
at each region to determine a vehicle routing policy. The 
objective of the control center is to detect an anomaly at 
any region in minimum time subject to a desired bound on 
the expected time between any two subsequent false alarms. 
Notice that the time required to detect an anomaly depends on 
the anomaly detection algorithm and the time vehicles take to 
travel the regions. Thus, the control center needs to minimize 
the anomaly detection time jointly over anomaly detection 
policies and vehicle routing policies. Our problem setup is 
shown in Fig. [T] 

We adopt the standard motion planning notation in ( [LaValle 



2006 1. We denote the k-th region by 7Zk, k G {1, . . . , n}, and 
r-th vehicle by Ur^r G {1, . . . , m}. Let the likelihood of an 
anomaly at region 1Zk be tt^ G (0, 1). We study the persistent 
surveillance problem under the following assumptions. 

Regarding the vehicles, we do not assume any specific 
dynamics and we assume that: 

(Al) each vehicle takes time dij to travel from region 1Zi to 

region Uj, ij G {l,...,n}; 
(A2) the sensors on each vehicle take a random time to 

collect an informative observatioij^from region 1Zkjk e 

^The vehicle routing policies designed in this paper also work for non- 
identical vehicles. We make this assumption for the convenience of analysis. 

^An informative observation may require the acquisition of several obser- 
vations from different locations at the same region. In this case the processing 
time equals the total time required to collect all these observations. 



{l,...,n}. 

Regarding the observations, we assume that: 

(A3) the observation collected by a vehicle from region IZk 
is sampled from probability density functions : R ^ 
R>o and : R ^ R>o, respectively, in the presence 
and in the absence of anomalies; 

(A4) for each k G {l,...,n}, probability density functions 
and /{? are non-identical with some non-zero proba- 
bility, and the two distributions have the same support; 

(A5) conditioned on the presence or absence of anomalies, the 
observations in each region are mutually independent; 
and 

(A6) observations in different regions are also mutually inde- 
pendent. 

Regarding the anomaly detection algorithm at the control 
center, we employ the cumulative sum (CUSUM) algorithm 
(see below) for anomaly detection at each region. In particular, 
we run n parallel CUSUM algorithms (one for each region) 
and declare an anomaly being present at a region as soon as 
a substantial evidence is present. We refer to such parallel 
CUSUM algorithms by ensemble CUSUM algorithm. 

Remark 1 (Knowledge of distributions): For the ease of 
presentation, we assume that the probability density func- 
tions in presence and absence of an anomaly are known. In 
general, only the probability density function in absence of 
any anomaly may be known, or both the probability density 
functions may be unknown. In the first case, the CUSUM 
algorithm can be replaced by the weighted CUSUM algorithm 
or the Generalized Likelihood Ratio (GLR) algorithm ( |Bas- 
|seville and Nikiforov[ |1993| ), while in the second case, it can 
be replaced by the robust minimax quickest change detection 
algorithm ( [Unnikrishnan et aL 2011 ). The ideas presented in 
this paper extend to these cases in a straightforward way. A 
related example is in Section [Vl] □ 

Remark 2 (Independence of observations): For the ease of 
presentation, we assume that the observations collected from 
each region are independent conditioned on the presence and 
absence of anomalies. In general, the observations may be 
dependent and the dependence can be captured through an 
appropriate hidden Markov model. If the observations can be 
modeled as a hidden Markov model, then the CUSUM like 
algorithm in ( Chen and Willett| 2000) can be used instead of 
the standard CUSUM algorithm. The analysis presented in this 
paper holds in this case as well but in an asymptotic sense, i.e., 
in the limit when a large number of observations are needed 
for anomaly detection. 

We also assumed that the observations collected from dif- 
ferent regions are mutually independent. Although the ideas in 
this paper also work when the observations at different regions 
are dependent, the performance can be improved with a slight 
modification in the procedure presented here (see Remark [4]). 
In this case the algorithm performance improves because each 
observation is now informative about more than one region. □ 

Regarding the vehicle routing policy, we propose the ran- 
domized routing policy, and the adaptive routing policy. In 
the randomized routing policy, each vehicle (i) selects a 
region from a stationary distribution, (ii) visits that region. 
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(iii) collects an evidence, and (iv) transmits this evidence 
to the control center and iterates this process endlessly. In 
the randomized routing policy, the evidence collected by the 
vehicles is not utilized to modify their routing policy. In other 
words, there is no feedback from the anomaly detection algo- 
rithm to the vehicle routing algorithm. In the adaptive routing 
policy, instead, the evidence collected by the vehicles is used 
to modify the routing policy, and thus, the loop between the 
vehicle routing algorithm and the anomaly detection algorithm 
is closed. The adaptive routing policy follows the same steps 
as in the randomized routing policy, with the exception that the 
distribution in step (i) is no longer stationary and is adapted 
based on the collected evidence. 

For brevity of notation, we will refer to the joint anomaly 
detection and vehicle routing policy comprising of the en- 
semble CUSUM algorithm and the randomized routing policy 
by randomized ensemble CUSUM algorithm. We will show 
that the randomized ensemble CUSUM algorithm provides a 
solution that is within a factor of optimality. Similarly, we 
refer to the joint anomaly detection and vehicle routing policy 
comprising of the ensemble CUSUM algorithm and adaptive 
routing policy by adaptive ensemble CUSUM algorithm. We 
will show that adaptive ensemble CUSUM algorithm makes 
the vehicles visit anomalous regions with high probability, 
and thus it improves upon the performance of the randomized 
ensemble CUSUM algorithm. The following standard defini- 
tion (Cover and Thomas] [1991 ) will be used in the remaining 
sections. 

Definition 1 (Kullback-Leibler divergence): Given two 
probability mass functions fi'-S^ M>o and f2 - S ^ ^>o, 
where S is some countable set, the Kullback-Leibler 
divergence V : x ^ RU{+oo} is defined by 



TABLE I 
List of symbols 



P(/l,/2)=% 



log 



fliX) 

f2{X) 



E 

a;esupp(/i) 



/i(a;)lo{ 



AM 

f2{xy 



where is the set of integrable functions, IE/J-] represents 
expected value with respect to /i, X is a random variable 
sampled from /i, and supp(/i) is the support of /i. □ 
It is known that (i) < /2) < +oo, (ii) the lower bound 

is achieved if and only if /i = /2 almost everywhere, and (iii) 
the upper bound is achieved if and only if the support of /2 is 
a strict subset of the support of /i. Observe that Assumption 
(A4) on the observations is equivalent to > 0, for 

each /c G {1, . . . , n}. 

We now introduce some notations that will be used through- 
out the paper. We denote the probability simplex in 
by A^_i, and the space of vehicle routing policies by Vt. 
For the processing time T/e, we let Tk denote its expected 
value. Consider m realizations of the processing time T^, 
we denote the expected value of the minimum of these m 
reaUzed values by T^-^^ist. Note that f^i-^^^ist = Tfe. We also 
define Tmax = max{f'/e | k e {l,...,n}} and Tmin = 
min{T/e | k G {l,...,n}}. We denote the Kullback-Leibler 
divergence between the probability density functions and 
by Vk. Finally, P^ax = max{P/e | A: G {1, . . . , n}} and 
Pmin = mm{Vk I /c G {1, . . . For the convenience of 

the reader, we have enlisted the notation in Table [H 



n 




number of regions 


m 




number of robots 


Uk^ke {1,. 


.,n} 


A;-th region 


Ur,r G {1,.. 




r-th vehicle 


TTA;, /e e {1, . . 


.,n} 


prior probability of anomaly at 1Zk 








dij,i,j G {1, . 


..,n} 


travel time between TZi and 1Zj 


Tk,k e {1,.. 


.,n} 


processing time at 7Zk 


fk,ke {1,.. 


■,n} 


expected processing time at 7Zk 


^max 




max{Tyt 1 k G {l,...,n}} 


^min 




min{Tfc | k G {l,...,n}} 


T~'m-smlst 




E[min{T^^\ . . . , T^'^^}], where 

Tj^ \ ' are m reahzations of Tk 


2^m-smlst 
min 




^in|j^m-smlst \ ke n}} 






set of sets of m arbitrary regions 


Tone 




min{E[min{t^, . . . , 1^}] C ^ where 
ti s are the processmg times at regions m ^ 


4°,fcs{i,.. 


.,n} 


pdf in absence of anomaly at 7lk 


fl,ke{i,.. 


.,n} 


pdf in presence of anomaly SLt TZk 


Vk,k e {1,. 


.,n} 


K-L divergence between and /° 


^max 




max{X>/i. 1 k G {1, . . . , n}} 






min{Vk 1 G {1, . . • ,n}} 






space of vehicle routing poHcies 


Nk : ^ ^ NU{+oo} 


observations required for detection at 1Zk 


Sk : ^ ^ M>o U{+oo} 


detection delay at 7Zk 


5avg : ^ M>o 


U{+cx)} 


T^k^i^kHSkiuj)] 


Supper : ^ — ^ ^>0 


U{+cx)} 


upper bound to 5avg 


rm-min 
^k 




inf{E[Sk(uj)] \uj en} 


rm-min 
"avg 




inf{5avg(^) 1 ^ G n} 






CUSUM statistic at TZj at r-th iteration 


7] 




CUSUM threshold 


fj 




e-v +77-1 


An-1 




probability simplex in 


q e An- 


1 


single vehicle randomized routing policy 




-1 


optimal q 




-1 


efficient q 


qm e A- 


-1 


m vehicle randomized routing policy 


^part ^ 


-1 


q^ with region partitioning 


a e An- 


1 


single vehicle adaptive routing policy 


apart e a;[ 


-1 


m vehicle adaptive routing poHcy 

with region partitioning 



Remark 3 (Randomized routing policy): The randomized 
routing policy samples regions to visit from a stationary 
distribution; this assumes that each region can be visited from 
another region in a single hop. While this is the case for 
aerial vehicles, it may not be true for ground vehicles. In 
the latter case, the motion from one region to another can be 
modeled as a Markov chain. The transition probabilities of this 
Markov chain can be designed to achieve a desired stationary 
distribution. This can optimally be done, for instance, by 



picking the fastest mixing Markov chain proposed in ( [Boyd 
et al.[ |2004| ) or heuristically by using the standard Metropolis- 
Hastings algorithm (jWasserman) |2004[). Related examples are 
presented in Section Vl] and |VII| It should be noted that 



5 



under the randomized routing poHcy, the desired stationary 
distribution of the Markov chain is fixed, and the Markov 
chain converges to this distribution exponentially. Thus, the 
policy designed using Markov chain is arbitrarily close to the 
desired policy. However, in the case of adaptive routing policy, 
the desired stationary distribution keeps on changing, and the 
performance of the Markov chain based policy depends on rate 
of convergence of the Markov chain and the rate of change of 
desired stationary distribution. □ 

III. Spatial Quickest Detection 

In this section we propose the ensemble CUSUM algorithm 
for the simultaneous quickest detection of anomalies in spa- 
tially distributed regions. We start by recalling the standard 
quickest change detection problem. Then we describe and 
characterize the ensemble CUSUM algorithm. 

A. Quickest change detection 

Consider a set of observations {^1,^2, • • where, for some 
z/, the observations {yi^ . . . ^y^-i} are i.i.d. with probabil- 
ity density function and {?/^, t/^^+i, . . .} are i.i.d. with 
probability density function The objective of the quickest 
change detection is to detect the change in the underlying 
distribution in minimum number of observations subject to 
a desired lower bound on the number of samples between 
two false alarms. Let N > u the observation at which 
the change is detected. The non-Bayesian quickest detection 
problem ( Poor and Hadjiliadis 2008 Siegmund 1985 ), is 
posed as 



minimize supEjy[A^ — v ^1\N > v] 
subject to Ejo[iV] > I/7, 



(1) 



where represents expected value with respect to the 

observations distribution at iteration u and 7 > is a small 
constant called /a/^^ alarm rate. 

An algorithmic solution to the minimization problem ([TJ 
is the cumulative sum (CUSUM) algorithm ( [Poor and Had-| 
jiliadis 2008| ), in which, at each iteration r G N, (i) an 
observation is collected, (ii) the statistic A^- = (At-_i + 



log 



with Aq = is computed, and (iii) a change 



is declared if A^- > r]. For a given threshold r], the false 
alarm rate and the worst expected number of observations for 
CUSUM algorithm are 



Ejo(7V) 



1 



1 



(2) 



The approximations in equation ^ are referred to as the 
Wald's approximations (Siegmund, 1985 ), and are known to be 
accurate for large values of the threshold r]. In the following, 
we assume that the chosen threshold is large enough and 
the expressions in equation ^ are exact. Let > be the 
uniform time duration between two iterations of the CUSUM 
algorithm. The expected detection delay 8, i.e., the expected 
time required to detect an anomaly after its appearance, 
satisfies E/i[^] = vEfi{N). 



Algorithm 1: Ensemble CUSUM Algorithm 

Input : threshold 77, pdfs , /c G {1, . . . , n} ; 

Output : decision on presence of an anomaly ; 

1 at time r receive observation for region TZ^; 

2 update the CUSUM statistic at each region: 



Ai 



k; 

ifj6{l, 



.,n}\{k}; 



3 if AJ > ?7 then change detected at region IZk '■> 

4 else wait for next observations and iterate. 



B. Ensemble CUSUM algorithm 

We run n parallel CUSUM algorithms (one for each region), 
and update the CUSUM statistic for region 1Zk only if an 
observation is received from region TZk- We refer to such 
parallel CUSUM algorithms by ensemble CUSUM algorithm 
(Algorithm [T]). Notice that an iteration of this algorithm is 
initiated by the collection of an observation. 

We are particularly interested in the performance of the 
ensemble CUSUM algorithm when the observations are col- 
lected by autonomous vehicles. In this case, the performance 
of the ensemble CUSUM algorithm is a function of the vehicle 
routing policy. For the ensemble CUSUM algorithm with 
autonomous vehicles collecting observation, let the number 
of iterations (collection of observations) required to detect 
an anomaly at region 1Zk be : ^ NU{+oo}, and 
let the detection delay, i.e., the time required to detect an 
anomaly, at region 1Zk he Sk : ^ ^ M>o U{+oo}, for each 
k e {l,...,n}, where ft is the space of vehicle routing 
policies. We also define average detection delay as follows: 

Definition 2 (Average detection delay): For any vector of 
weights {wi^ . . . ^Wn) G A^-i, define the average detection 
delay (^avg - ^ ^ ^>o U{+oo} for the ensemble CUSUM 
algorithm with autonomous vehicles collecting observations 
by 



(3) 



k=i 



For the ensemble CUSUM algorithm with m vehicles 
collecting observation, define ^^-"^^^ and S'^'^^^ by 



rm-min 



= inf{E[4(a;)] \ujen}, and 
= inf{(5avg(^) I Co' G ^1}, 



respectively. Note that S'^'^^'^ and S'^'^^'^ are lower bounds for 
the expected detection delay and average detection delay at 
region 1Zk, respectively, independently of the routing policy. 
Let f] = -\- 1] — 1. We now state lower bounds on 
the performance of the ensemble CUSUM algorithm with 
autonomous vehicles collecting observations. 

Lemma 3 (Global lower bound): The following statements 
hold for the ensemble CUSUM algorithm with m vehicles 
collecting information: 



6 



(i) the lower bound S^'^^"^ for the expected detection delay 
at region 1Zk satisfies 



> 



mVk 



(ii) the lower bound ^^g"^^^ for the average detection delay 



satisfies 



> 



m-smlst 
min 



mVrr 



where T^^^^'' = min{Tf -^^^^^ | A: G {1, . . . , n}}. 

Proof: We start by establishing the first statement. We 
note that a lower bound on the expected detection delay 
at region 1Zk is obtained if all the vehicles always stay at 
region 1Zk- Since, each observation is collected from region 

the number of iterations of the ensemble CUSUM al- 
gorithm required to detect an anomaly at region 1Zk satisfies 
E[Nk] = fj/Vk- Let Tl(h) be realized value of the processing 
time of vehicle Ur at its 6-th observation. It follows that 
2^m-smist,^^) = min{T^(6) | r G {l,...,m}} is a lower 
bound on the processing time of each vehicle for its 6-th 
observation. Further, jj^-smistj^^j identically distributed for 
each h and ^[T^-'^^^'^h)] 



j^m-smist Consider a modified 



stochastic process where the realized processing time of each 
vehicle for its 6th observation in j^-smistj^^j Indeed, such a 
stochastic process underestimates the time required to collect 
each observation and, hence, provides a lower bound to 
the expected detection delay. Therefore, the detection delay 
satisfies the following bound 

\Nk/m] 

Skioj) > Yl T^"'^^'\b), for each uj e Q. 

6=1 



It follows from Wald's identity ( |Resnick[ |1999 |), that 

This proves the first statement. 

The second statement follows from Definition [2] and the first 
statement. ■ 

Remark 4 (Dependence across regions): We assumed that 
the observations collected from different regions are mutually 
independent. If the observations from different regions are 
dependent, then, at each iteration, instead of updating only one 
CUSUM statistic, the CUSUM statistic at each region should 
be updated with an appropriate marginal distribution. □ 

IV. Randomized Ensemble CUSUM Algorithm 

We now study the persistent surveillance problem under 
randomized ensemble CUSUM algorithm. First, we derive 
an exact expression for the expected detection delay for 
the randomized ensemble CUSUM algorithm with a single 
vehicle, and use the derived expressions to develop an ef- 
ficient stationary policy for a single vehicle. Second, we 
develop a lower bound on the expected detection delay for 
the randomized ensemble CUSUM algorithm with multiple 
vehicles, and develop a generic partitioning policy that (i) 
constructs a complete and disjoint m-partition of the regions, 
(ii) allocates one partition each to a vehicle, and (iii) lets each 



vehicle survey its assigned region with some single vehicle 
policy. Finally, we show that the partitioning policy where 
each vehicle implements the efficient stationary policy in its 
regions is within a factor of an optimal policy. 

A. Analysis for single vehicle 

Consider the randomized ensemble CUSUM algorithm with 
a single vehicle. Let qk G [0, 1] denote the probability to select 
region 1Zk^ and let = (g'l, . . . , qn) G A^_i. Let the threshold 
for the CUSUM algorithm at each region be uniform and equal 
to 7^ > 0. We note that for the randomized ensemble CUSUM 
algorithm with a single vehicle the space of vehicle routing 
policies is = A^-i. 

Theorem 4 (Single vehicle randomized ensemble CUSUM): 
For the randomized ensemble CUSUM algorithm with a 
single vehicle and stationary routing policy q G A^-i, the 
following statements hold: 

(i) the number of observations Nk{q) required to detect a 
change at region 1Zk satisfies 



(ii) the detection delay Sk{q) at region 1Zk satisfies 

n n n 

% [4(q)] = ( E + E E mdi^Yfl \Nk{q)]. 

i=l i=l j=l 

Proof: Let r G {1, . . . , A/'^} be the iterations at which 
the vehicle collects and sends information about the regions, 
where Nk denotes the iteration at which an anomaly is 
detected at region 1Zk- Let the log likelihood ratio at region 
at iteration r be A^. We have 

Xk ^ fllf^, with probability qk, 

[O, with probability I — qk- 

Therefore, conditioned on the presence of an anomaly, 
{A^}t-^i^ are i.i.d., and 

Ef^[X'^]=qkVk. 

The remaining proof of the first statement follows similar to 



the proof for CUSUM in ( |Siegmund| [1985). 

To prove the second statement, note that the information 
aggregation time T^^^ comprises of the processing time and 
the travel time. At an iteration the vehicle is at region 1Zi 
with probability qi and picks region 1Zj with probability qj. 
Additionally, the vehicle travels between the two regions in dij 
units of time. Thus, the average travel time at each iteration 
is 

n n 

E[Ttravel] = ^^QiQjdjj- 

i=l j=l 

Hence, the expected information aggregation time at each 
iteration is 

n n n 
E[T^Sr] = E[r,ravel + ^process] = E E + E 

i=l j=l i=\ 



7 




Fig. 2. Level- sets of the objective function in problem J5j. It can be seen 
that the level sets are not convex. 



Let {Tr^^}^^|i TVfc}. be the information aggregation times at 
each iteration. We have that 6k = X^^i Tr^\ and it follows 
from Wald's identity ( |Resnick| |1999| ), that 



This completes the proof of the statement. ■ 

B. Design for single vehicle 

Our objective is to design a stationary policy that si- 
multaneously minimizes the detection delay at each region, 
that is, to design a stationary policy that minimizes each 
term in {Si{q)^ . . . ^Sn{q)) simultaneously. For this multiple- 
objective optimization problem, we construct a single aggre- 
gate objective function as the average detection delay. After 
incorporating the expressions for the expected detection delays 
derived in Theorem [4j the average detection delay becomes 

Hj Th Th Th 

^ i=l i=l j=l 

where Wk = 'T^k/ (Yl^^i^i) is the weight on the expected 
detection delay at region IZk and is the prior probability 
of an anomaly being present at region IZk- Our objective is to 
solve the average detection delay minimization problem: 



minimize (^avg(9)- 



(5) 



In general, the objective function (^avg is non-convex. For 
instance, let n = 3, and consider the level sets of (5avg on the 
two dimensional probability simplex (Fig. [2]). It can be seen 
that the level sets are non-convex, yet there exists a unique 
critical point and it corresponds to a minimum. We now state 
the following conjecture about the average detection delay: 

Conjecture 5 (Single vehicle optimal stationary policy): 
For the randomized ensemble CUSUM algorithm with a 
single vehicle, the average detection delay function ^avg has 
a unique critical point at which the minimum of (5avg is 
achieved. □. 

In the Appendix we provide probabilistic guarantees that, 
for a particular stochastic model of the parameters in (5avg, 
with at least confidence level 99.99% and probability at least 
99%, the optimization problem ([5]) has a unique critical point 
at which the minimum is achieved. Such a minimum can be 



computed via standard gradient-descent methods (|Boyd and 



Vandenberghel 2004) 



We now construct an upper bound for the expected detection 
delay. We will show that minimization of this upper bound 
yields a policy that is within a factor of an optimal policy. 
From equation (|4]), we define the upper bound ^, 



>o U{+oo} as 



upper 



n _ 

^avg(<7) < (^upper(<7) = ( ) (?max + C^max), 

where (imax = max{(ii^- | i, j G {1, . . . , n}}. 

Theorem 6 (Single vehicle efficient stationary policy): 
The following statements hold for the randomized ensemble 
CUSUM algorithm with single vehicle: 

(i) the upper bound on the expected detection delay satisfies 



mm 



-'upper 



n 

(^) = (E 



H~ ^max ) 7 



k=l 



and the minimum is achieved at defined by 



ke {l,...,n}; 



(ii) the average detection delay satisfies the following lower 
bound ^ 

k=l ^ ^ 

for all q G A^-i; 

(iii) the stationary policy q^ is within a factor of optimal, 
that is 



^avg 



■ d„ 



Trr 



< n- 



and 



where gr* is an optimal stationary policy; 
(iv) the expected detection delay at region IZk under policy 



satisfy 



rl-min 
^k 



WkVn 



Proof: We start by establishing the first statement. It 
follows from the stationarity conditions on the Lagrangian that 
the minimizer gr^ of (^upper satisfy ql cx ^/wkf]/T)k, for each 
k G {l,...,n}. Incorporating this fact into YI^^iqI = 1 
yields the expression for q^. The expression for (Supper (9^) can 
be verified by substituting the expression for q^ into (Supper- 

To prove the second statement, we construct a lower bound 
Slower • ^ ^>o U{+oo} to the average detection delay 

^avg defined by (5iower(<7) = Y.k=i^kf]fmin/T>kqk- It can be 
verified that ^lower also achieves its minimum at q^ , and 

Slower (g^"^) = {Yl \lW) ^ 
k^l V 

We note that 

^low&r{q^) < ^lower(Q*) < ^avg(Q*) < ^avg(g'),Vgr G A^-l- 
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Thus, the second statement follows. 

To prove the first part of the third statement, we note that 

^upper(9^)- 

Therefore, the policy is within (Supper 
(^max + c^max)/^min factor of Optimal Stationary poHcy. 
To prove the second part of the third statement, we note 



Algorithm 2: Partitioning Algorithm 



max V-*- max 



V ■ T ■ 

-^mm-^ mm 



< n- 



7 V 

'f'max -^max 



-mm -^mm 

where the last inequality follows from the fact: max{y^ 
. . . + \ wi ^ . . . ^ Wn = 1} = \fn. 

To establish the last statement, we note that 

max; 



s:l-min 



< 



< 



< 



WkVn 



/Wr, 



(^max H~ ^max) / Til^i 



This concludes the proof of the theorem. ■ 
In the following, we would refer to as the single vehicle 
efficient stationary policy. 

Remark 5 (Efficient stationary policy): As opposed to the 
average detection delay S^yg, the upper bound (^upper does not 
depend upon any travel time information. Then, our efficient 
policy does not take this information into account, and it may 
not be optimal. Instead, an optimal policy allocates higher 
visiting probabilities to regions that are located more centrally 
in the environment. We resort to the efficient policy because 
(i) if the problem ([Sj does not admit a unique minimum, then 
the optimal policy can not be computed efficiently; and (ii) 
the efficient policy has an intuitive, tractable, and closed form 
expression. □ 

C. Analysis for multiple vehicles 

We now consider the randomized ensemble CUSUM with 
m > 1 vehicles. In this setting the vehicles operate in 
an asynchronous fashion. This asynchronicity, which did not 
occur in the single vehicle case, is due to (i) different travel 
times between two different pair of regions, and (ii) different 
realized value of processing time at each iteration. Such an 
asynchronous operation makes the time durations between two 
subsequent iterations non-identically distributed and makes it 
difficult to obtain closed form expressions for the expected 
detection delay at each region. 

Motivated by the above discussion, we determine a lower 
bound on the expected detection delay for the randomized 
ensemble CUSUM algorithm with multiple vehicles. Let q^ = 
{Qi^ ' ' ' ^Qn) ^ denote the stationary policy for vehicle 

Ur, i.e., the vector of probabilities of selecting different regions 
for vehicle Ur, and let = (g^, . . . , g^) G A5^_^. We note 
that for the randomized ensemble CUSUM algorithm with m 
vehicles the space of vehicle routing policies is = A5^_i. We 



Input : vehicles {Ui , . . . , Um}, regions 7Z = {7Zi , . . . , IZn}, 
a single vehicle routing policy; 

Require : n > m ; 

Output : a m-partition of the regions ; 

1 partition TZ into m arbitrary subsets {S^}re{i,...,m} 

with cardinalities rir < \n/m],r G {!,..., m} ; 

2 allocate vehicle Ur to subset S^, for each r G {!,..., m}; 

3 implement the single vehicle efficient stationary policy in each subset. 



construct a lower bound on the processing times at different 
regions for different vehicles in the following way. Let S be the 
set of all the sets with cardinality m in which each entry is an 
arbitrarily chosen region; equivalently, S = {T^i, . . . ,7^n}^- 
Let a realization of the processing times at the regions in a set 
^ G S be 4 7 • • • 7 • define a lower bound Tone to the 

expected value of the minimum of the processing times at m 
arbitrary regions as Tone = min{E[min{t^, • . • , t^}] | ^ G S}. 

Theorem 7 (Multi-vehicle randomized ensemble CUSUM): 
For the randomized ensemble CUSUM algorithm with 
m vehicles and stationary region selection policies q^, 
r G {1, . . . , m}, the detection delay at region 1Zk satisfies: 



Proof: We construct a modified stochastic process to 
determine a lower bound on the expected detection delay. For 
the randomized ensemble CUSUM algorithm with multiple 
vehicles, let be the the processing time for the vehicle 
Ur during its 6-th visit to any region. We assume that the 
sampling time for each vehicle at its 6-th visit in the modified 
process is min{tj, . . . , t^}. Therefore, the sampling time for 
the modified process is the same at each region. Further, 
it is identically distributed for each visit and has expected 
value greater than or equal to Tone- We further assume that 
the distances between the regions are zero. Such a process 
underestimates the processing and travel time required to 
collect each observation in the randomized ensemble CUSUM 
algorithm. Hence, the expected detection delay for this process 
provides a lower bound to the expected detection delay for 
randomized ensemble CUSUM algorithm. Further, for this 
process the vehicles operate synchronously and the expected 
value of the likelihood ratio at region k at each iteration is 
YlT^i Qk^ifk^ fk)- remainder of the proof follows similar 
to the proof for single vehicle case in Theorem |4] ■ 

D. Design for multiple vehicles 

We now design an efficient stationary policy for random- 
ized ensemble CUSUM algorithm with multiple vehicles. We 
propose an algorithm that partitions the set of regions into m 
subsets, allocates one vehicle to each subset, and implements 
our single vehicle efficient stationary policy in each subset. 
This procedure is formally defined in Algorithm [2] 

Let the subset of regions allocated to vehicle Z//^ be 5^, r G 
{l,...,m}. We will denote the elements of subset 5^ by 
Sl,i G {1, . . . , n^}. Let q^^^^ G A^_i be a stationary routing 
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policy under the partitioning algorithm that implements single 
vehicle efficient stationary policy in each partition. We define 
the weights in equation ([3]) by = tt^/ Trj, where 
ttI is the prior probability of anomaly at region 1Zk- Let 
^min = min{^i, . . . , Wn} and ^max = max{^i, . . . , Wn}. 
We now analyze the performance of the partitioning algorithm 
and show that it is within a factor of optimal. 

Theorem 8 (Performance of the partitioning policy): For 
the partitioning algorithm with m vehicles and n regions that 
implements the single vehicle efficient stationary policy in 
each partition, the following statements hold: 

(i) the average detection delay under partitioning policy 
satisfies the following upper bound 



^avg(^part) < ^ 



max \ '-f'max 



(ii) the average detection delay satisfies the following lower 
bound 



n I r 

W9J>(EV^)^ 



k=l 



for any G A! 



n-l' 



(iii) the stationary policy q^^^^ is within a factor of optimal, 
and 



^avg(^Jart) ^ ^Wma.^ {T^ 



^avg (^rn) 



^avg(^part) . 2 

— - — -. — < m 



.)Vn 



-, and 



^m-smlst 



Vn 



where q*^ is optimal stationary policy; 
(iv) the expected detection delay at region IZk under the 
stationary policy q^^^^^ satisfies 



part 



n 



rm-mm 



< 



T-im-smlst 



Proof: We start by establishing the first statement. We 
note that under the partitioning policy, the maximum number 
of regions a vehicle serves is [n/m] . It follows from Theo- 
rem [6] that for vehicle Ur and the associated partition 5^, the 
average detection delay is upper bounded by 



^avg(<7part) < ( a/^) ^(^max + ^^ma 
i=l ^ * 
r on -1 2 ono/i <; (^max H~ <^max ) 



:) 



the overall average detection delay satisfies 
avg(5jart) ^ ^^avg(9part)- ^his establishes the first statement. 



Therefore, 

To prove the second statement, we utilize the lower bounds 
obtained in Theorem |7] and construct a lower bound to the 
average detection delay S^^^^^ : A5^_^ ^R>oU{+oo} defined 

by <wer(^m) = ELi (^fcfone/Er=i • ^an be Verified 
that 

Vk / m 



mm 



"^lower 



(Qm) 



Algorithm 3: Single Vehicle Adaptive Ensemble CUSUM 



Input 
Output 

1 set Aj = 0, for all j E {1, 
while true do 

2 



parameters 77, V^, pdfs f^, f^, for each /c G {1, 
decision on anomaly at each region ; 

n}, and r = 1; 



set new prior tt^ — -l ^A, 



set Qk : 



e^-/(l + e^T), for each k G {1, 
, for each /c G {1, . . . , n}; 



sample a region from probability distribution (qi 

collect sample i/r from region k; 

update the CUSUM statistic at each region 



Ai 



(At 



if j = k- 

if jG{l,... 



,n}\W; 



if Aj > r; then 

anomaly detected at region TZj^; 
set Aj = 0; 
set r = r + 1 ; 



We now establish the first part of the third statement. Note 
that 



^avg(^Jart) ^ \n/m\^ ^max {T^ 



.)Vn 



Vn 



< 



.)Vn 



Wrr 



Vn 



k=l 



where the last inequality follows from the fact that 
{\n/m])/{n/m) < 2. 

The remainder of the proof follows similar to the proof of 
Theorem [6l ■ 

V. Adaptive ensemble CUSUM Algorithm 

The stationary vehicle routing policy does not utilize the 
real-time information regarding the likelihood of anomalies 
at the regions. We now develop an adaptive policy that 
incorporates the anomaly likelihood information provided by 
the anomaly detection algorithm. We consider the CUSUM 
statistic at a region as a measure of the likelihood of an 
anomaly at that region, and utilize it at each iteration to 
design new prior probability of an anomaly for each region. At 
each iteration, we adapt the efficient stationary policy using 
this new prior probability. This procedure results in higher 
probability of visiting an anomalous region and, consequently, 
it improves the performance of our efficient stationary policy. 
In Section |Vl| we provide numerical evidence showing that 
the adaptive ensemble CUSUM algorithm improves the per- 
formance of randomized ensemble CUSUM algorithm. 

Our adaptive ensemble CUSUM algorithm is formally pre- 
sented in Algorithm [3] for the single vehicle case. For the case 
of multiple vehicles we resort to the partitioning Algorithm [2] 
that implements the single vehicle adaptive ensemble CUSUM 
Algorithm [3] in each partition. Let us denote the adaptive 
routing policy for a single vehicle by a and the policy obtained 
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from the partitioning algorithm that implements single vehicle 
adaptive routing policy in each partition by apart- We now 
analyze the performance of the adaptive ensemble CUSUM 
algorithm. Since, the probability to visit any region varies with 
time in the adaptive ensemble CUSUM algorithm, we need to 
determine the number of iterations between two consecutive 
visit to a region, i.e., the number of iterations for the recurrence 
of the region. We first derive a bound on the expected number 
of samples to be drawn from a time- varying probability vector 
for the recurrence of a particular state. 

Lemma 9 (Mean observations for region recurrence): 
Consider a sequence where is sampled from a 

probability vector G A^-i. If the k\h entry of satisfy 
Pk ^ {^k^h), for each r G N and some ak^h ^ (0? 1). 
then the number of iterations 1^ for the recurrence of state k 
satisfy E[4] < Pu/al 

Proof: The terms of the sequence are statisti- 

cally independent. Further, the probability mass function is 
arbitrary. Therefore, the bound on the expected iterations for 
the first occurrence of state k is also a bound on the subsequent 
recurrence of state k. The expected number of iterations for 
first occurrence of region k are 



\{{^-p^k)<PkY.'i^i-akr- 



Pk/Oil- 



This establishes the statement. ■ 
We utilize this upper bound on the expected number of 
iterations for recurrence of a region to derive performance 
metrics for the adaptive ensemble CUSUM algorithm. We 
now derive an upper bound on the expected detection delay 
at each region for adaptive ensemble CUSUM algorithm. We 
derive these bounds for the expected evolution of the CUSUM 
statistic at each region. 

Theorem 10 (Adaptive ensemble CUSUM algorithm): 
Consider the expected evolution of the CUSUM statistic at 
each region. For the partitioning algorithm that implements 
single vehicle adaptive ensemble CUSUM algorithm 
(Algorithm [3]) in each subset of the partition, the following 
statement holds: 



E[4(apart)] < 

, {\n/m\ - 



2{\n/m\ - l)e^/2yp^(i _ 



V^min 



^min(l - e- 



<)• 



Proof: We start by deriving expression for a single 
vehicle. Let the number of iterations between the (j — l)th 
and jth visit to region IZk be . 

Let the observation during the jth visit to region IZk be 
yj and the CUSUM statistic at region IZk after the visit be 
. It follows that the probability to visit region IZk between 
(j — l)th and jth visit is greater than 



Vk 



Therefore, it follows from Lemma [9] that 

ni-] < (1 + (n - l)e^''-'^U)/^^V^, 



Note that C>y = max{0, Cj=_i + log(/i(%)//fe°(t/i))}. Since, 
maximum of two convex function is a convex function, it 
follows from Jensen inequality ( |Resnick[ [1999), that 

Therefore, E[Cj] > jV^ and for expected evolution of the 
CUSUM statistics 

ni^] < (1 + (n - l)e^^-^^--i)^^)/V^i^/^min)'. 

Therefore, the total number of iterations A^^ required to 
collect N'^^ observations at region TZ}. satisfy 



nNk{a)\Nt'] = ^(1 + (n - l)e(''-0-i)^'')/ V^'fc/I'min)' 
2{n-l)e^'^^k{l-e-^'^<''^) 



= Nt"- -r ; 

' v^(l-e-^^/2) 

(n-l)^e''I>fc(l-e-P^^*^) 

Note that the number of observations AT*** required at region 
TZk satisfy E[A^jg''*] = fj/Vk-lt follows from Jensen's inequal- 
ity that 



E[Nk{a)] < 



Jl_ 2{n - l)e''/^v^(l - e-^/^) 
(1 - e-^^) 
(n-l )^e^Pfc(l - e-^) 



'mm(l-e-^0 ■ 
Since the expected time required to collect each evidence is 
smaller T^ax + <^max, it follows that 

E[4(a)] < (Tmax+^imax)E[7Vfe(a)]. 

The expression for the partitioning policy that implement sin- 
gle vehicle adaptive routing policy in each partition follow by 
substituting [n/m] in the above expressions. This completes 
the proof of the theorem. ■ 
Remark 6 (Performance bound): The bound derived in 



Theorem 10 is very conservative. Indeed, it assumes the 
CUSUM statistic at each region to be fixed at its maximum 
value T], except for the region in consideration. This is practi- 
cally never the case. In fact, if at some iteration the CUSUM 
statistic is close to r], then it is highly likely that the vehicle 
visits that region at the next iteration, so that the updated 
statistic crosses the threshold r] and resets to zero. □ 

VI. Numerical Results 

We now elucidate on the concepts developed in this paper 
through some numerical examples. We first validate the ex- 



pressions for expected detection delay obtained in Section IV 



Example 1 (Expected detection delay): Consider a set of 
4 regions surveyed by a single vehicle. Let the loca- 
tion of the regions be (10, 0), (5, 0), (0, 5), and (0,10), re- 
spectively. The vector of processing times at each region 
is (1,2,3,4). Under the nominal conditions, the observa- 
tions at each region are sampled from normal distribu- 
tions A/'(0,l),A/'(0, 1.33), A/'(0, 1.67) and A/'(0,2), respec- 
tively, while under anomalous conditions, the observations 



11 




Threshold 



(b) Average detection delay 

Fig. 3. Expected and average detection delay. Solid black line with 
dots and black x, respectively, represent the theoretical expression and the 
value obtained by Monte-Carlo simulations under stationary policy q = 
[ 0.2 0.25 0.25 0.3 ]. Dashed green line and green triangles, respectively, 
represent the theoretical expression and the value obtained by Monte-Carlo 
simulations under stationary policy q = [ 0.5 0.2 0.2 0.1 ]. Solid 
red line and red diamonds, respectively, represent the theoretical expression 
and the value obtained by Monte-Carlo simulations under stationary policy 
q = [0.85 0.05 0.05 0.05]. 

are sampled from normal distributions with unit mean and 
same variance as in nominal case. Let the prior probability of 
anomaly at each region be 0.5. An anomaly appears at each 
region at time 50, 200, 350, and 500, respectively. Assuming 
that the vehicle is holonomic and moves at unit speed, the 
expected detection delay at region IZi and the average de- 
tection delay are shown in Fig. [3] It can be seen that the 
theoretical expressions provide a lower bound to the expected 
detection delay obtained through Monte-Carlo simulations. 
This phenomenon is attributed to the Wald's approximation. □ 
We remarked in Section [ll] that if each region cannot be 
reached from another region in a single hop, then a fastest 
mixing Markov chain (FMMC) with the desired stationary 
distribution can be constructed. Consider a set of regions 
modeled by the graph Q = {V^£), where V is the set of nodes 
(each node corresponds to a region) and £ is the set of edges 
representing the connectivity of the regions. The transition 
matrix of the FMMC P G R^^^ with a desired stationary 
distribution q G A^_i can be determined by solving the 
following convex minimization problem ( (Boyd et slI] |2004|: 

minimize WQ'^^PQ'^^ - Qroot^rLlb 
subject to PI = 1 

QP = P^Q 

Pij > 0, for each (z, j) G £ 
Pij = 0, for each (i, j) ^ f , 

where Q is a diagonal matrix with diagonal q, q^^^^ = 
{\/Q1j ' ' ' j -s/Qm)^ and 1 is the vector of all ones. We now 
demonstrate the effectiveness of FMMC in our setup. 

Example 2 (Effectiveness of FMMC): Consider the same 
set of data as in Example [T] We study the expected and 
average detection delay for randomized ensemble CUSUM 
algorithm when the regions to visit are sampled from the 




Threshold 



(a) Expected detection delay at region TZ\ 



^ i2oor 




Threshold 

(b) Average detection delay 

Fig. 4. Expected and average detection delay for uniform stationary policy. 
The solid black line represents the theoretical expression. The black x, red 
diamonds, and green triangles, respectively, represent the values obtained by 
Monte-Carlo simulations for all-to-all, line, and ring connection topology. For 
the line and ring topologies, the region to visit at each iteration is sampled 
from the fastest mixing Markov chain with the desired stationary distribution. 




Fig. 5. Average detection delay for a single vehicle. The solid red line, 
the dashed green line, and the solid black line with dots represent efficient, 
uniform, and optimal stationary policies, respectively. 

FMMC. The expected and average detection delay for all- 
to-all connection topology, line connection topology and ring 
connection topology are shown in Fig. |4] It can be seen that 
the performance under all three topologies is remarkably close 
to each other. □ 

We now study the performance of the (numerically com- 
puted) optimal and our efficient stationary policies for the 
single vehicle randomized ensemble CUSUM algorithm. 

Example 3 (Single vehicle optimal stationary policy): 
For the same set of data as in Example [T] we now study 
the performance of the uniform, the (numerically computed) 
optimal and our efficient stationary routing policies. A 
comparison is shown in Fig. [5] Notice that the performance 
of the optimal and efficient stationary policy is extremely 
close to each other. □ 

We now study the performance of the optimal, partitioning 
and uniform stationary policies for randomized ensemble 
CUSUM algorithm with multiple vehicles. 

Example 4 (Multiple-vehicle optimal stationary policy): 
Consider a set of 6 regions surveyed by 3 vehicles. Let the 
regions be located at (10, 0), (5, 0), (0, 5), (0, 10), (0, 0) and 
(5,5). Let the processing time at each region be unitary. 
Under nominal conditions, the observations at each region 
are sampled from normal distributions A/'(0, 1), A/'(0, 1.4), 
A/'(0,1.8), A/'(0,2.2), A/'(0,2.6) and A/'(0,3), respectively. 
Under anomalous conditions, the observations are sampled 
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▲ 




Fig. 6. Average detection delay for 3 vehicles surveying 6 regions. The 
green triangles represent the policy in which each vehicle surveys each region 
uniformly. The red diamonds and black x represent the partitioning policy in 
which each vehicle implements the single vehicle efficient stationary policy 
and the single vehicle optimal stationary policy, respectively. 
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from normal distributions with unit mean and same variance 
as in the nominal case. Let the prior probability of anomaly 
at each region be 0.5. An anomaly appears at each region 
at time 25,35,45,55,65 and 75, respectively. Assuming that 
the vehicles are holonomic and moves at unitary speed, the 
average detection delay for the uniform stationary policy for 
each vehicle, the partitioning policy in which each vehicle 
implements single vehicle efficient stationary policy in each 
subset of the partition, and the partitioning policy in which 
each vehicle implements single vehicle optimal stationary 
policy in each subset of the partition is shown in Fig. [6] □ 

We now study the performance of the adaptive ensemble 
CUSUM algorithm, and we numerically show that it improves 
the performance of our stationary policy. 

Example 5 (Adaptive ensemble CUSUM algorithm): 
Consider the same set of regions as in Example [T] Let the 
processing time at each region be unitary. The observations at 
each region are sampled from normal distributions A/'(0,cr^) 
and A/'(l, cr^), in nominal and anomalous conditions, 
respectively. Under the nominal conditions at each region 
and cr^ = 1, a sample evolution of the adaptive ensemble 
CUSUM algorithm is shown in Fig. [7ja). The anomaly 
appears at regions 7^3, and IZ4 at time 100, 300, 
and 500, respectively. Under these anomalous conditions 
and cr^ = 1, a sample evolution of the adaptive ensemble 
CUSUM algorithm is shown in Fig. [TJb). It can be seen that 
the adaptive ensemble algorithm samples a region with high 
likelihood of anomaly with high probability, and, hence, it 
improves upon the performance of the stationary policy. 

We now study the expected detection delay under adaptive 
ensemble CUSUM algorithm and compare it with the efficient 
stationary policy. The anomaly at each region appears at time 
50, 200, 350 and 500, respectively. The expected detection 
delay obtained by Monte-Carlo simulations for = 1 and 
different thresholds is shown in Fig. [SJa). It can be seen that 
the adaptive policy improves the detection delay significantly 
over the efficient stationary policy for large thresholds. It 
should be noted that the detection delay minimization is most 
needed at large thresholds because the detection delay is 
already low at small thresholds. Furthermore, frequent false 
alarms are encountered at low thresholds and hence, low 
thresholds are not typically chosen. The expected detection 
delay obtained by Monte-Carlo simulations for different value 
of cr^ and threshold 7^ = 5 is shown in Fig.jsj^b). Note that for 
a given value of a^, the Kullback-Leibler divergence between 



(a) CUSUM statistic and vehicle routing probabilities under nominal conditions 
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(b) CUSUM statistic and vehicle routing probabilities under anomalous conditions 

Fig. 7. Sample evolution of the adaptive ensemble CUSUM algorithm. The 
dashed-dotted blue line, dashed green line, solid red line and solid black 
line with dots represent data from regions IZi , 7^2, '^s and 7^4, respectively. 
The solid brown horizontal line represents the threshold. The vehicle routing 
probability is a function of the likelihood of anomaly at each region. As 
the likelihood of an anomaly being present at a region increases, also the 
probability to survey that region increases. Anomalies appear at region 7l2, 
IZ3 and 7^4 at times 100, 300 and 500, respectively. Once an anomaly is 
detected, it is removed and the statistic is reset to zero. 



A/'(l,cr^) and A/'(0,cr^) is l/2cr^. It can be seen that the 
adaptive policy improves the performance of the stationary 
policy for each value of noise. □ 
We now apply the adaptive ensemble CUSUM algorithm 
to a more general scenario where the anomalous distribution 
is not completely known. As remarked in Section in this 
case, the CUSUM algorithm should be replaced with the GLR 
algorithm. Given the nominal probability density function 
and the anomalous probability density function fl{-\0) 
parameterized by 6> G 6 C M^, for some ^ G N, the GLR 
algorithm ("Basse ville and Nikiforov 1993 ), works identically 
to the CUSUM algorithm, except that the CUSUM statistic is 
replaced by the statistic 



max supVlog%#. 



Example 6 (Generalized Likelihood Ratio): For the same 
set of data as in Example [5] assume that there are three types 
of potential anomalies at each region. Since any combination 
of these anomalies can occur simultaneously, there are 7 
potential distributions under anomalous conditions. We char- 
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Threshold 

(a) Expected detection delay as a function of threshold 




(K-L Divergence)" 

(b) Expected detection delay as a function of KL divergence 

Fig. 8. Performance of the adaptive ensemble CUSUM algorithm. The solid 
black line represents the theoretical expected detection delay for the efficient 
stationary policy and the black x represent the expected detection delay for 
the adaptive ensemble CUSUM algorithm. 



acterize these distributions as different hypothesis and assume 
that the observations under each hypothesis h G {1, . . . , 8} 
are sampled from a normal distribution with mean /i/^ and 
covariances S^. Let 
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We picked region 1Zi as non-anomalous, while hypothesis 4, 
6, and 8 were true at regions 7^2, '^3, and respectively. 
The Kullback-Leibler divergence at a region was chosen as the 
minimum of all possible Kullback-Leibler divergences at that 
region. A sample evolution of the adaptive ensemble CUSUM 
algorithm with GLR statistic replacing the CUSUM statistic is 
shown in Fig[9fa). It can be seen the performance is similar to 
the performance in Example |5] As an additional ramification of 
this algorithm, we also get the likelihood of each hypothesis at 
each region. It can be seen in Fig[9jb) that the true hypothesis 
at each region corresponds to the hypothesis with maximum 
likelihood. □ 

VII. Experimental Results 

We first detail our implementation of the algorithms us- 
ing the Player/Stage robot control software package and the 
specifics of our robot hardware. We then present the results of 
the experiment. 

Robot hardware: We use Erratic mobile robots from Videre 
Design shown in Fig. [TO] The robot platform has a roughly 
square footprint (40cm x 37cm), with two differential drive 
wheels and a single rear caster. Each robot carries an on- 
board computer with a L8Ghz Core 2 Duo processor, 1 GB of 
memory, and 802. llg wireless communication. For navigation 
and localization, each robot is equipped with a Hokuyo URG- 
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(a) GLR statistic under anomalous conditions 



Hypothesis 

(b) Normalized likelihood of each hypothesis 

Fig. 9. Sample evolution of the adaptive ensemble CUSUM algorithm with 
GLR statistic. The dashed-dotted blue line, dashed green line, solid red line 
and solid black line with dots represent data from regions IZi, IZ2, ^3 and 
7Z4, respectively. The solid brown horizontal line represents the threshold. The 
vehicle routing probability is a function of the likelihood of anomaly at each 
region. As the likelihood of an anomaly being present at a region increases, 
also the probability to survey that region increases. Anomalies appear at region 
^2, ^3 and IZ4 at times 100, 300 and 500, respectively. Once an anomaly 
is detected, it is removed and the statistic is reset to zero. The true hypothesis 
at each region corresponds to the hypothesis with maximum likelihood 
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Fig. 10. Erratic mobile robot with URG-04LX laser rangefinder. 



04LX laser rangefinder. The rangefinder scans 683 points over 
240° at lOHz with a range of 5.6 meters. 

Localization: We use the amcl driver in Player which 
implements Adaptive Monte-Carlo Localization ( [Thrun et"aLl 
^2001). The physical robots are provided with a map of our 
lab with a 15cm resolution and told their starting pose within 



the map (Fig. [TT]). We set an initial pose standard deviation of 
0.9m in position and 12° in orientation, and request updated 
localization based on 50 of the sensors range measurements 
for each change of 2cm in robot position or 2° in orientation. 
We use the most likely pose estimate by amcl as the location 
of the robot. 

Navigation: Each robot uses the snd driver in Player for 
the Smooth Nearness Diagram navigation ( [Durham and BuUo] 
2008). For the hardware, we set the robot radius parameter 
to 22cm, obstacle avoidance distance to 0.5m, and maximum 
speed to 0.2m/s. We let a robot achieve its target when it is 
within 10cm of the target. 
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Fig. 11. This figure shows a map of our lab together with our surveillance 
configuration. Three erratic robots survey the selected 8 regions (black dots), 
which have been partitioned among the robots. Regions 1, 2, and 3 are also 
considered in Fig. [12] where we report the statistics of our detection algorithm. 



Experiment setup: For our experiment we employed our 
team of 3 Erratic robots to survey our laboratory. As in 
Fig. [TT] a set of 8 important regions have been chosen and 
partitioned among the robots. Each robot surveys its assigned 
regions. In particular, each robot implements the single robot 
adaptive ensemble CUSUM algorithm in its regions. Notice 
that Robot 1 cannot travel from region 1 to region 3 in a 
single hop. Therefore, Robot 1 selects the regions according 
to a Markov chain with desired stationary distribution. This 
Markov chain was constructed using the Metropolis-Hastings 
algorithm. In particular, for a set of regions modeled as a graph 
Q = (VjS), to achieve a desired stationary routing policy q, 
the Metropolis-Hastings algorithm ( ,Wasserman| |2Q04| ), picks 
the transition matrix P with entries: 



0, 



mm 



En 
k= 



Qidj J 



if (i, j) G E and i 7^ j, 
if (i, j) G f and z = j, 



where di is the number of regions that can be visited from 
region IZi. 

Observations (in the form of pictures) are collected by a 
robot each time a region is visited. In order to have a more 
realistic experiment, we map each location in our lab to a 
region in our campus. Then, each time a robot visit a region in 
our lab, a picture of a certain region in our campus is selected 
as observation (see Fig. [12]). Pictures have been collected prior 
to the experiment. 

Finally, in order to demonstrate the effectiveness of our 
anomaly detection algorithm, some pictures from regions 2 
and 3 have been manually modified to contain an anomalous 
pattern; see Fig. [13] Anomalous pictures are collected by 



Robot 1 at some pre-specified time instants (the detection 
algorithm, however, does not make use of this information). 

Probability density function estimation: In order to imple- 
ment our adaptive ensemble CUSUM algorithm, the proba- 
bility density functions of the observations at the regions in 
presence and absence of an anomaly need to be estimated. 




Fig. 12. A snapshot of our surveillance experiment, where three robots 
survey six locations in our lab (Fig. [TT). In this fi sure we show the three 
regions assigned to the first robot. Each region correspond to a part of 
our campus, and observations are taken accordingly. Notice that Region 3 
contains an anomaly (black smoke), and that the CUSUM statistics, which 
are updated upon collection of observations, reveal the anomaly (green peak). 
The transition probabilities are updated according to our adaptive ensemble 
CUSUM algorithm. 
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Region 2 with anomaly 




Region 3 



Region 3 with anomaly 




Fig. 13. This figure shows sample pictures from Region 2 and Region 3, 
both with and without the anomaly to be detected. 



For this task, we first collect sample images, and we register 



them in order to align their coordinates ( [Radke et al.| |2QQ5| ). 
We then select a reference image, and compute the difference 
between the sample pictures and the reference image. Then, 
we obtain a coarse representation of each difference image 
by dividing the image into blocks. For each difference image, 
we create a vector containing the mean value of each block, 
and we compute the mean and standard deviation of these 
vectors. Finally, we fit a normal distribution to represent the 
collected nominal data. In order to obtain a probability density 
distribution of the images with anomalies, we manually modify 
the nominal images, and we repeat the same procedure as in 
the nominal case. 

Experiment results: The results of our experiment are 



illustrated in Fig. 12 Fig. 13 and in the multimedia extension 
available at Ihttp : / / www . i j rr . org, From the CUSUM 
statistics we note that the anomalies in Region 2 and Region 
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3 are both detected: indeed both the red curve and the 
green curve pass the decision threshold. We also note that 
few observations are necessary to detect the anomaly. Since 
the robots successfully survey the given environment despite 
sensor and modeling uncertainties due to real hardware, we 
conclude that our modeling assumptions in Section [ll] are not 
restrictive. 

VIII. Conclusions 

In this paper we studied a spatial quickest detection prob- 
lem in which multiple vehicles surveil a set of regions to 
detect anomalies in minimum time. We developed a novel 
ensemble CUSUM algorithm to detect an anomaly in any of 
the regions. A stochastic vehicle routing policy was adopted 
in which the vehicle samples the next region to visit from 
a probability vector. In particular, we studied (i) stationary 
policy: the probability vector is a constant function of time; 
and (ii) adaptive policy: the probability vector is adapted with 
time based on the collected observations. We designed an 
efficient stationary policy that depends on the travel time of the 
vehicles, the processing time required to collect information 
at each region, and the anomaly detection difficulty at each 
region. In adaptive policy, we modified the efficient stationary 
policy at each iteration to ensure that the regions with high 
likelihood of anomaly are visited with high probability, and 
thus, improved upon the performance of the stationary policy. 
We also mentioned the methods that extend the ideas in this 
paper immediately to the scenario in which the distributions of 
the observations in presence and absence of anomaly are not 
completely known, but belong to some parametrized family, 
or to the scenario in which the observations collected from 
each region are not independent (e.g., in the case of dynamic 
anomalies). 

There are several possible extensions of the ideas considered 
here. First, in the case of dependent observations at each 
region, the current method assumes known distributions in 
presence and absence of anomalies. An interesting direction 
is to design quickest detection strategies that are robust to 
the uncertainties in these distributions. Second, the anomalies 
considered in this paper are always contained in the same 
region. It would be of interest to consider anomalies that can 
move from one region to another. Third, the policy presented 
in this paper considers an arbitrary partition that satisfy some 
cardinality constraints. It is of interest to come up with smarter 
partitioning policies that take into consideration the travel 
times, and the difficulty of detection at each region. Last, to 
construct the fastest mixing Markov chain with desired sta- 
tionary distribution, we relied on time-homogeneous Markov 
chains. A time varying Markov chain may achieve a faster 
convergence to the desired stationary distribution ( [Grace and] 
Baillieul 20051 ). This is also an interesting direction to be 
pursued. 

Appendix 

A. Probabilistic guarantee to the uniqueness of critical point 

We now provide probabilistic guarantee for Conjecture |5] 
The average detection delay for a single vehicle under a 



stationary policy q is 



^avg 



n n 



Qj dij ^ , 



where Vi = Wifj/Vi for each i G {l,...,n}.A local minimum 
of ^avg can be can be found by substituting = 1 — Yl^Zi Qj ' 
and then running the gradient descent algorithm from some 
initial point Qq G A^_i on the resulting objective function. 

Let V = (^1, . . . , Vn) and T = (Ti, . . . , T^). We assume 
that the parameters {v^T^ D^n} in a given instance of op- 
timization problem ^ and the chosen initial point Qq are 
realizations of random variables sampled from some space 
JC. For a given realization G /C, let the realized value of 
the parameters be {v{k,)^T{k,)^ D{K,)^n{K,)}, and the chosen 
initial point be g^o(^)- The associated optimization problem is: 



minimize S^yg{q \k,)^ 



(A-1) 



where, for a given realization k, e JC, (^avg(' | ^) • ^n{K)-i 
]R>o U{+oo} is defined by 



^avg(<7 I ^) = ( XI "^^) (IZ Qi^M^Y^ qiQjdijif^fj 

i=i i=i i=i j=i 

For a given realization hz, let gd{- \ n) : A^(^)_i ^n{K)-i 
be the function that determines the outcome of the gradi- 
ent descent algorithm applied to the function obtained by 
substituting q^^^) = 1 - ^Ij in (^avg(<7 | a^). In other 

words, the gradient descent algorithm starting from point go(^) 
converges to the point gd{qo{hz) \ hz). Consider A^i realizations 

{A^i,...,/^ArJ e JC^\ Let gr^ptiniai(^) = ^^(^ln(/.) I 
and define 

7 = max{||^d((7o(/^s)|/^s)-goptimai(^^)ll I {l,...,Ni}}. 

It is known ( |Calafiore cTH] [20TT] ) that if Ni > 
— (log for some G ]0, 1[, then, with at least 

confidence 1 — i^i, it holds 

n(/^)— 1 I \\9d{qo {hz) I hz) 9optimal(^)ll <7}) 

for any realization hz e JC. 

We sample the following quantities: the value n as uni- 
formly distributed in {3,.. .,12}; each coordinate of the n 
regions in two dimensional space from the normal distribu- 
tion with mean and variance 100; the value T^, for each 
i G {l,...,n}, from the half normal distribution with mean 
and variance 100; and the value Vi, for each i G {l,...,n}, 
uniformly from ]0, 1[. For a realized value of n, we chose 
qQ uniformly in A^_i. Let the matrix D be the Euclidean 
distance matrix between the n sampled regions. 

We considered Ni = 1000 realizations of the parameters 
{v^T^ D^n} and initial value qQ. The sample sizes were 
determined for /ii = 0.01 and ui = 10 The value of 
7 obtained was 10~^. Consequently, the gradient descent 
algorithm for the optimization problem ^ starting from any 
feasible point yields the same solution with high probability. 
In other words, with at least confidence level 99.99% and 
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probability at least 99%, the optimization problem ^ has a 
unique critical point at which the minimum is achieved. 
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